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Derivation of F ormulations 1 and 1 A 

of Farassat 

F. Farassat 

NASA Langley Research Center, Hampton, Virginia 

Summary 

Formulations 1 and 1A are the solutions of the Ffowcs Williams-Hawkings (FW-H) equation with 
surface sources only when the surface moves at subsonic speed. Both formulations have been success- 
fully used for helicopter rotor and propeller noise prediction for many years although we now recom- 
mend using Formulation 1A for this purpose. Formulation 1 has an observer time derivative that is taken 
numerically, and thus, increasing execution time on a computer and reducing the accuracy of the results. 
After some discussion of the Green’s function of the wave equation, we derive Formulation 1 which is 
the basis of deriving Formulation 1A. We will then show how to take this observer time derivative 
analytically to get Formulation 1A. We give here the most detailed derivation of these formulations. 
Once you see the whole derivation, you will ask yourself why you did not do it yourself! 

1- Primary References 

Formulation 1 was first published in: 

1- F, Farassat: Theory of Noise Generation From Moving Bodies With an Application to Helicop- 
ter Rotors, NASA Technical Report R-451, 1975 

However, the derivation of this formulation in this reference is very difficult. Later, I found a much 
simpler derivation which we will follow here. It was published in: 

2- F. Farassat: Linear acoustic formulas for calculation of rotating blade noise, AIAA Journal, 
19(9), 1981, 1122-1130 

Formulation 1A was first published in: 

3- F. Farassat, G. P. Succi, A review of propeller discrete frequency noise prediction technology 
with emphasis on two current methods for time domain calculations, 1980, Journal of Sound and 
Vibration , 71(3), 399-419 
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In deriving the formulation in the above reference, I assumed, as it was common then, that the blade 
surface did not move out of the disc plane. This assumption was removed by Kenneth S. Brentner and 
was published in the following excellent NASA Technical Memorandum: 

4- Kenneth S, Brentner, Prediction of Helicopter Discrete Frequency Rotor Noise- A Computer 
Program Incorporating Realistic Blade Motions and Advanced Formulation, October 1986, 
NASA TM 87721 

We refer to the result published in this paper as Formulation 1A. At the time of writing Reference 3, 
we at Langley had not started numbering and identifying our formulations. In hindsight, it was a good 
idea to identify each of our formulations because different codes use different results. Our numbering 
system has worked for identifying what formulation is used in a code. Formulation 1A has been used in 
helicopter rotor noise prediction codes WOPWOP and PSU-WOPWOP, and together with our super- 
sonic Formulation 3, in our advanced propeller noise prediction code ASSPIN (Advanced Subsonic 
and Supersonic Propeller Induced Noise). 

You will need much mathematical maturity to follow what we discuss here. To understand the mathemat- 
ics behind our work, I strongly recommend reading Reference 1, above, and the following two NASA 
publications: 

5- F. Farassat, Introduction to Generalized Functions With Applications in Aerodynamics and 
Aeroacoustics, May 1994 (Corrected April 1996), NASA Technical Paper 3428 

6- F. Farassat: The Kirchhoff Formulas for Moving Surfaces in Aeroacoustics - The Subsonic and 
Supersonic Cases, NASA Technical Memorandum 110285, September 1996 

Reference 6 was written to clarify and fill in some gaps in the mathematical analysis of Reference 5. It 
should be read together with the latter reference. 

2- The Ffowcs Williams-Hawkings (FW-H) 
Equation 


This equation was first published in: 

7- J. E. Ffowcs Williams and D. L. Hawkings: Sound generated by turbulence and surfaces in 
arbitrary motion, Philosophical Transactions of the Royal Society, A264, 1969, 321-342 

This paper is very difficult to read because of the high level of mathematics needed to follow the 
authors’ reasoning. There are many original ideas in this paper such as the idea of imbedding a problem 
in a larger domain to use the available Green’s function of the larger domain to solve the original prob- 
lem. Another idea put forward by Ffowcs Williams and Hawkings is that conservation laws in differen- 
tial form are also valid when all ordinary derivatives are viewed as generalized derivatives. This has 
important implications about the jump conditions across flow discontinuities. I have given the back- 
ground mathematics needed to understand this paper in References 1, 5 and 6. 

Ffowcs Williams-Hawkings (FW-H) Equation as originally proposed in Reference 7 above: 
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□ 2 P' = -“[po v„ 6(f)] 
ot 


d 

dxi 


[ P rii 6(f)] + 


d 2 

d Xj dxj 


[ H(f ) Ty] 


( 1 ) 


Here and elsewhere in the paper, we will use the summation convention on the repeated index. In this 
equation D 2 is the wave or D’Alembertian operator in three dimensional space. The moving surface is 
described by /( x, t) = 0 such that vf = n, n is the unit outward normal. This assumption implies that 
/ > 0 outside the moving surface (see Figure 1). Also p' - c 2 p' - c 2 (p - pf ) , c and po are the speed 
of sound and density in the undisturbed medium, respectively. Note that p' can only be interpreted as 
the acoustic pressure if p' / Po « 1 . The symbols v n , p and 7)y = p Uj Uj - cr ,-y + (// - c 2 p 1 ) d (/ are the 
local normal velocity of the surface, the local gage pressure on the surface (in fact p - po), and the 
Lighthill stress tensor, respectively. In the definition of the Lighthill stress tensor, <Xy is the viscous 
stress tensor and <5,y is the Kronecker delta. The Heaviside and the Dirac delta functions are denoted 
H{f) and 6(f ) , respectively. In the second term on the right of eq. (1), we have neglected the viscous 
shear force over the blade surface acting on the fluid. 



f(xj)< 0 


Figure 1- The definition of the moving surface implicitly as f(x, t) = 0. Note that v / = n where n 
is the unit outward normal to the surface. 

We note that we have artificially converted a nonlinear problem of noise generation by a moving surface 
to a linear problem by using the acoustic analogy. All the nonlinearities are lumped into the Lighthill 
stress tensor which is assumed known from near field aerodynamic calculations. When we started work- 
ing on helicopter and propeller noise in the early seventies, because of the limitations of digital comput- 
ers, the most we could expect from aerodynamic calculations was the blade surface pressure. For this 
reason, using some physical reasoning, we neglected the quadrupole volume sources in FW-H equation 
and concentrated on development of formulations for the prediction of thickness and loading noise. 
Later on, as computers became more powerful, we included quadrupoles in our noise prediction. There 
was, however, another theoretical advance which led to the use of purely surface sources to which we 
will turn next. 

It was Ffowcs Williams himself who proposed to use a penetrable (porous or permeable) data surface to 
account for nonlinearities in the vicinity of a moving surface. We again assume that the penetrable 
surface defined by f (x, t) = 0 and the fluid velocity is denoted by u. The FW-H equation for penetra- 
ble (permeable, porous) data surface, FW-H p( j s , is: 
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n c\ c\2 

a 2 c 2 p' = □ V = — [p„ C/„] <S(/) - — ft 6(f)] + — — [ T,j H(f)] 

Ot OXi OXjOXj 

We have used the following notations in the above equation: 


U n 



V„ + 


P Un 
Po 


( 2 ) 

(3) 


Li = p 6 ij rij + p Uj (u n - v„) 


(4) 


where <5y is the Kronecker delta. As in the case of eq. (1), in the first term on the right of eq. (4), we 
have neglected the viscous shear force over the data surface acting on the fluid exterior to the surface. 
The philosophy behind using FW-H pc js is to locate the data surface / = 0 to enclose a moving surface, 
in such a way that all quadrupoles producing non-negligible noise are included within this surface. 
Therefore, no volume integration of the quadrupoles outside the data surface is necessary. High resolu- 
tion CFD calculation is performed in the near field region (including turbulence simulation if broadband 
noise prediction is required). The data surface used for the acoustic calculation should be located within 
the region of high resolution CFD computation. The optimal location of the data surface, particularly 
when vortices cross the surface, is still the subject of research. One would like this surface to be as small 
as possible, because of the computer intensive nature of aerodynamic and turbulence simulation that 
require fine grid sizes and small time steps. 

We will be concerned with the solution of the following two wave equations: 


Thickness Noise Equation (5) 


Loading Noise Equation (6) 


n 2 p' T = — [p 0 v„ 6(f)] 

□ 2 P' L = --~[pni8{f)] 


The source terms on the right of eqs. (5) and (6) are known also as monopole and dipole sources, respec- 
tively. This terminology is misleading for sources on a moving surface because the radiation patterns 
calculated from the above equations are not similar to those from stationary monopoles and dipoles. In 
fact, we can show mathematically that the thickness noise is equivalent to the loading noise from a 
uniform surface pressure distribution of magnitude po c 1 . For this reason we avoid using the terminol- 
ogy of monopole and dipole sources here. 

Remark 1- We will now give an example of a moving surface defined implicitly by / = 0. A sphere of 
radius a moving at the speed v along the jq-axis is described by the equation: 

fi = (X| - v t ) 2 + x\ + x\ - a 2 =0 (7) 

Note that f 1 > 0 outside the surface, and f) < 0 inside it. Now let us test the gradient of this function. 
We see that 

X\ - vt X 2 X 2 \ 


v/i = (2 (xi - vt), 2x 2 , 2 x 2 )+ n = 


a 


a 


a 


( 8 ) 
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We have \V f[ \ = 2 ^ (x\ - v t) 2 + x 2 + x 2 . Therefore, we redefine the moving sphere by the equation: 
fix, t) = ■ = — a/ (X| - V t) 2 + 4 + 4 - 


Iv/il 


yfixi-vt) 


2 +x\+x 2 


(9) 


We give the following general rule: if V/ =£ n on the surface f - 0, then redefine the same surface by 
the implicit equation //|v/| = 0. We can then show that v(//|v/|) = n on this surface. The proof is 
simple. 

Note that a surface can be defined implicitly in more than one way satisfying the condition v / = n. For 
example, for the above surface a much simpler representation in implicit form satisfying this condition 

is f{x, t) - ^ (x\ - v t) 2 + x\ + x 2 - a = 0 . 

We have not said anything about why we require the condition |V/| = 1 on the surface. The FW-H 
equation as it was originally written in Reference 7 had the term |V/| in both the thickness and loading 
terms. We have been assuming |v/| = 1 in writing the FW-H equation in all of our publications for 
many years. To begin with, we notice that the term | v f\ does not appear in both Formulations 1 and 1A. 
In the process of the derivation of some of our more advanced formulations where the source time and 
space derivatives of |V/| were required, I noticed that all terms associated with these derivatives can- 
celed out exactly in the final result. After considerable search for the reason or reasons behind this 
cancellation, I discovered that |v/| = 1 could be assumed to be true on any surface described implicitly 
from the beginning of the derivation of FW-H equation. This explained the reason for cancellation of the 
terms associated with the derivatives of |V/| in our advanced formulations. Although this assumption is 
of little value in the derivation of Formulation 1, it reduces the algebra somewhat in the derivation of 
Formulation 1A and enormously in our advanced formulations. (End of Remark 1) 


3- Background Material and the Details of 
the Derivation 


3.1- Green’s Function of the Wave Equation in Unbounded 
Three Dimensional Space 

The Green’s function of the wave equation in the unbounded three dimensional space is: 


G(x, t\ y, r) = 


0 T > t 

6(T-t + r/c)/4nr r^t 


( 10 ) 


where r = \ x - y \ . Here (x, t) and (y, r) are the observer and the source space-time variables, respec- 
tively. In the above equation, the symbol 6( ■ ) stands for the Dirac delta function which is the most well- 
known generalized function. We usually use the following symbol: g = T — t + r/c. We will discuss in 
Remark 3 below how we visualize and interpret the function g = 0 . See also References 1 and 5 above. 

The Green’s function given by eq. (10) is also known as the free-space Green’s function. 
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Remark 2- It is a good idea to remember the following two simple results because they are frequently 
encountered in algebraic manipulations: 


dr dr 

= r: and = - r,- 

dx{ dy t 


( 11 ) 


where r, is the component of the unit radiation vector (x - y)/r. (End of Remark 2) 

Remark 3- The surface g = T- t + r/c = 0 can be visualized as follows. Keep the observer space-time 
variables fixed. Then let us first rewrite the equation of this surface as | x - y | = c (t - r) . In the source 
space-time variables, this surface is simply a sphere with center at x and radius equal to c (t - t). Note 
that the Green’s function of the wave equation is nonzero when r t, so that as the source time 
increases from r - -oo to r = t, the radius of this sphere shrinks from infinitely large value to zero. For 
this reason this sphere is known as the collapsing sphere. The rate of contraction of the radius of this 
sphere is the speed of sound c . 

Mathematically, the collapsing sphere is the characteristic cone of the wave equation with the vertex at 
(x, t) and the cone pointing towards the past. Causality rules out using the part of the cone pointing 
towards the future. Can you think of a reason why this surface is called a cone? See Reference 6, above 
for the answer. (End of Remark 3) 


3.2- Solution of the Wave Equation With Sources on a Mov- 
ing Surface 

We are interested in the solution of the equation: 


□ V= Q(x,t)S(f) 


( 12 ) 


We will derive the solution of this equation using the free-space Green’s function given by eq. (10). 
This will show the interconnection between the variables that confuses novices in the field of aeroacous - 
tics. The formal solution of eq. (12) is: 

4np'(x,t) = J^Q(y, r)S(f) dydr (13) 


The limits of the integral in this equation are given below: 

r*t r*00 rsOO rsOQ 

dydr= I I dydr= I I I I dy\dy 2 dyi,dT (14) 

— OO'-'IR 1 — oo — oo — oo G — oo 

We emphasize at this point that the x - frame and the y - frame are fixed to the undisturbed medium. 
For the next analytic step, such a frame is not the most appropriate frame to interpret the solution of eq. 
(12) as we will show below. It is important to remember that from now on, the observer space-time 
variables (x, t) are kept fixed in all algebraic manipulations. Therefore, for all practical purposes, we are 
dealing with four variables (y, r) . 

For problems of interest in aeroacoustics, such as propeller and helicopter rotor noise prediction, one 
can always describe the surface (generally a blade) in a frame fixed relative to the surface. We will call 
this frame the rj - frame . Such a frame is used by the manufacturer to describe the design of the blade 
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surface to the technicians who build the blade. We call the variable ij the Lagrangian variable of a 
point on the moving surface. If we mark a point on the surface, say by a red dot, then we have essen- 
tially identified a point with the fixed variable // on the moving surface. In general, in the //-frame the 
moving surface is time independent although this does not have to be so in the following mathematical 
manipulations. The thing to realize at this point is that once the motion of the blade is specified, then the 
trajectory of a point on the surface described by a fixed ij is specified in space-time in the y- frame, i.e., 
in space. This trajectory is given by the relation: 


y = y(n, r) 


(15) 


The inverse transformation is given by: 

r/ = r/(y, r) 


(16) 


Note that the equation of the moving surface f(y, r) = 0 in the //-frame is f(y (ij, r), r) = /(//, r) = 0. In 
practice, / = 0 is independent of time, i.e., the surface is described as / (//), and this is what we assume 
here. Later on we will see that it is customary to use / when we should properly use f for designating 
the moving surface in an analytic expression (e.g., see eq. (29)). This is an abuse of notation which is 
acceptable for the reason that we will discuss in the paragraph following eq. (31). You should keep in 
mind that, in general, the analytic expressions for / and / are different. For example, in the example of 

a moving (rigid) sphere in Remark 1, we have f(y, t) = \J (>’i - v t) 2 + y\ + y\ - a = 0 while 

/ 67) = V^i + di + 73 - a = 0 . Here, we are assuming that the //-frame has its origin at the center of the 
sphere with its axes parallel to corresponding axes of the y-frame. However, this distinction does not 
matter to us here because we are not able to integrate any of the acoustic integrals here in closed form 
for any nontrivial moving surface. We always evaluate our surface integrals numerically by finite differ- 
ence scheme after we divide a data or a blade surface into panels. 

For problems of interest to us in aeroacoustics, the transformations described by eqs. (15) and (16) are 
isometric , i.e., distance preserving, because they involve translations and rotations only. For this reason, 
the Jacobians of transformation are unity, that is, we have 


detf^l 

= 1 

and 

detf^Ul 

V dr/ ) 



\dy) 


(17) 


If you have trouble accepting this, then use the translation transformation y = // + v r for a constant 
velocity v used in Remark 1 (where v = v e\ and e\ is the basis vector along t}\ - axis) to test the valid- 
ity of eq. (17). We now go back to the integration of the Dirac delta functions of eq. (13). We first use 
the transformation y -» // in this equation to get 

4 7T p' (x, t) = 

f Q(y(* h T), t) 6(J) ------- - - drj dr = f f Q (r/, r) S(f) dr/ d t (18) 

J r\det(dr//dy)\ J-ooJ/r 1 r 

where we have defined Q (//, r) = Q{y{ij , r), r). Next we use the transformation r -» g. The Jacobian of 
this transformation is dg / 8 t. Here the variable // is kept fixed in this partial differentiation. We will do 
the algebra in detail below. First we have 
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T ~t + 


x- yin, T) 

c 


The partial differentiation with respect to variable r gives 


d_g_ 

dr 


1 dr dyj r ; - v ; - 

c dyj dr c 


1 -M r 


(19) 


( 20 ) 


where M r =r jVj/ c is the Mach number of the point // in the radiation direction at the time r. Here r, is 
the component of unit radiation vector (x - y) / r and v ( - = dy i (//, r)/ <9r is the component of the velocity 
v of the point rj with respect to the y - frame fixed to the undisturbed medium. 


Using the above results in eq. (18) and assuming that s is a small positive number, we get 

4 np'(x,t)= f Q ( 1 1 , r) 6(f) - dj]dg = 

J r\dg/dr\ 



8(g) 


r II -Mr 


d gdrj 


-f 

JlR 3 


Q ( n , t) 

( r\l -M r \ 


8(f) 


dij 


8 = 0 


( 21 ) 


Note that the limits of the inside integral (with respect to g) of the expression after the second equality 
sign are from -s to e (e > 0) because d(g) could only contribute to the integral in this region. In what 
follows, we will drop the absolute value sign around |1 - M r \ because for surfaces moving subsonically, 
we always have 1 - M r > 0. The expression 1 - M r is known as the Doppler factor. Let us now inter- 
pret what exactly the condition g = 0 imposes on us. Remember that we have used the transformation 
r -> g. This means that, from eq. (19), the condition g - 0 makes the source time dependent on the other 
variables (x, t; rj) . This function is found analytically by solving for r from the equation 

g = t - t + fS-L — 11 _ o ((x, t) kept fixed) (22) 

c 


Before we go further, let us see what this source time signifies. Remember that when we fix tj, it means 
that we have marked a fixed point on the moving surface, say, with a red dot. Therefore, the trajectory 
of this red dot in space-time is known from the relation y = y(T],r). That is, given any source time r, 
we know where the red dot is on its trajectory in space. But from eq. (22), we have 

I x - y (tj, t) | = c (t - t) (23) 


Let us write r e = t(x, t\ //) and y e = y (//, r e ). We will shortly say what the subscript e stands for. With 
these notations, eq. (23) can be written as 

r e = \x-y e \ = c(t-T e ) (24) 

Clearly, the source time r e is the emission time, y e = y (//. ry ) is the emission position, and r e is the 
emission distance of the source point r/ to the observer position x. Both the emission position and the 
emission distance are viewed in the y - frame fixed to the undisturbed medium. See Figure 2 for an 
illustration of some of the terms we use here. It can be shown that for a subsonically moving source, 
there is only one emission time for a given source point r\ . In practical problems of rotating blade noise 
prediction, we cannot find a closed form analytic expression for r = T e {x, t\ tj) because eq. (23) is a 
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transcendental equation involving sines and cosines. However, we can find r e numerically easily by a 
shooting technique. Roughly speaking, here is what we do. At the observer time t, we know where the 
source (the red dot) position is and its trajectory in space. This position is called the visual position of 
the source. Walk in small time steps along the trajectory back in (source) time, i.e., for r < t, and each 
time check to see if eq. (23) is satisfied. You may be overshooting the emission position y(j], r e ). In 
that case, walk along the trajectory towards the visual position of the source in smaller time steps till 
you find the emission time and position satisfying eq. (23). In practice, this method can be refined 
considerably to speed up the computation of the emission time. 



Figure 2- Sketch of the trajectory of a source point rj (the red dot) as seen by an observer fixed to 
the medium ( x - or y-frame) and the definitions of visual and emission positions of the source 

With the above notations, eq. (21) can be written as: 


An p' (x, t) = 



QO h T e ) 
r e (\~M re ) 


$\f (.rf)]chj 


(25) 


where we have defined M r> _ = M(i /, r e ) ■ r e and M(i /, r e ) = v(i/, r e )/ c. Now note that 

Q(r], T e ) 
r e (\-M re ) 


= l/r (x, t\ 7 /) 


(26) 
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where the function if (x, t\ //) shows the dependence on variables of the expression on the left of eq. 
(26). We can show that for an arbitrary integrable function q(y), we have the following easily remem- 
bered and beautiful result (See References 5 and 6): 


(assuming that | V f\ = 1) (27) 


f q(y)8(f)dy= | 

q(y) dS 

Jk 3 j. 

f = o 


Note that y is a dummy variable here. Therefore, integrating the delta function in eq. (25) gives 


4 71 p' (x, t) = 



Q0 h T e ) 
r e (1 - M re ) 


dS 


(28) 


It is customary, but confusing and not entirely correct, to write the above result as: 



(29) 


where the subscript ret stands for the retarded time. The reason that this result is not entirely correct is 
that the expression on the left of eq. (26) is not obtained simply by replacing the source time r in the 
integrand of eq. (29) by t - r/c which is what the subscript ret implies. Although it is always true that 
r e - t - r e / c, r e itself is given by the expression | x - y(ij, T e ( x , t; //)) | which shows that the retarded 
time notation of eq. (29), and even the use of the source variable y in its integrand, is not permissible. 
Also note that in this equation we have used the notation / = 0 instead of the correct notation f = 0 . In 
summary, for a moving surface, we define: 


Q (y, t) 


Lr(l 


M r ) 


J ret 


Q (y, r) 

_ r (1 - M r ) 


(x; T], r) 


A T = T e (X, f,Tj) 


Q(J1, T e ) 

r e a~Mr e ) 


(30) 


where in the second square brackets, Q and M are functions of (//, r), and r and r are functions of 
(x; i], t) and then the whole thing is evaluated at the emission time. Other symbols are defined as fol- 
lows: 


Q0 h Te ) = Q(y(T], Te), T e ), r e = \x~y e \, = j(//, T e ), M, e = M{j], T e )-r e (31) 

This means that the left side of eq. (30) is a shorthand notation for its right side utilized to avoid the 
introduction of many new symbols, such as rj, f , Q(j], T e ), etc. However, as you will see below, when 
we seek the solution of FW-H equation, we always end up with integrals of the type in eq. (28). We will 
use the shorthand notation of eq. (29) repeatedly below without warning. Unfortunately, the notation of 
eq. (29) is the source of considerable confusion if one does not understand that the integrand of this 
equation is not to be interpreted by the conventional meaning of retarded time notation for stationary 
sources. For this reason, we have given here a detailed derivation of the solution of eq. (12). 

We emphasize here that, in general, we cannot get analytic expressions for the quantities in eq. (31) for 
most problems of interest in aeroacoustics. However, all these quantities can be found numerically with 
arbitrary precision. Combining the analytic solution of FW-H equation with the power of modern digital 
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computers allows us to use the exact geometry and kinematics of rotating blade machinery. This 
approach has been enormously successful in all areas of aeroacoustics. 

Remark 4- When a source moves rectilinearly at uniform speed, then the source time can be found in 
closed form by solving eq. (22). Furthermore, we discover that we have one emission time for subsonic 
motion and two emission times for supersonic motion of the source. (End of Remark 4) 

3.3- Derivation of Formulation 1 

We will now derive Formulation 1. As will be seen, for thickness noise, we will have no trouble writing 
the solution. Based on what we know about the solution of eq. (12), the solution of eq. (5), the thick- 
ness noise part of Formulation 1 can be written as: 


, d 

Anp T (x, t)= — 

f 

f= o 

Po v„ 

Lr(l - M r )\ 

dS 

ret 


We will say more about the mathematics behind keeping the time derivatives outside the integral sign 
later (See Remark 5). 

We will now derive the loading noise part of Formulation 1. We know that for an observer in the far 
field, we have the following approximation: 


_d_ ^ _ h_ d_ 
dxj c dt 


( 33 ) 


Can one derive the exact result? The affirmative answer was found by Farassat in Reference 1. We will 
give here the derivation in Reference 2 which is shorter and more elegant. Fet us start by writing the 
formal solution of eq. (6) using the free-space Green’s function (See Remark 6): 

4 7t p' L (x, t) = J' p rii 6(f) d y d r = - prii 6(f) ^ ^dydr ( 34 ) 


Now we will use the following identity which can be proven by carrying out the differentiation on both 
sides: 


(> ( <%) \ _ _ 1 d_( r j 6(g) \ _ r 6(g) 

dxi 1 r j c dt\ r ) r 2 

If we take the partial derivatives on both sides, we get 


d 1 

( <%) ) 

1 M'(g) M(g) 1 d | 

f r t S(g) ) 

j 1 M'(g) 

dxi 

(~i~ ) 

1 

! ^ 
1 ^ 

<N 

s*. 

! 

i ^ 

1 

( r J 

' c r 


( 35 ) 


( 36 ) 


This proves the identity of eq. (35). Fet us now use this identity on the right of eq. (34) and then take the 
observer time derivative out of the first integral. We obtain: 
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4 7T p' L (X,t)=-^j'p n t 6(f) 


d_ ( riS(g) 

dt 1 r 


dydr + 


J' P n i ?i 


8(f) 


8(g) 


d y d t = 


i or „ 

I p n; r ; 

c dt J y 


8(f) dydT + J priori 6(f) dydr 


(37) 


See Remark 6 about an important mathematical point in relation to this equation. Again using the 
solution of eq. (12), we write the solution of eq. (6) which is the loading noise part of Formulation 1 
as: 


4 7T p' L (x, t) = 


1 d 
c dt 


f 

J/=o 


p cos 6 

L r (T Z ~M r ) \ 


dS + 

f f 

ret •J 

/= o L 


p cos 6 

= 0 L r 2 (1 - M r )\ 


dS 


ret 


(38) 


where cos 0 = tij r l , i.e., 0 is the local angle between normal to the surface and radiation direction r at 
the emission time. 

Formulation 1 of Farassat is the sum of eqs. (32) and (38): 


47 T p' (x, t) = 4 7t (p' T (x, t) + p' L (x, t)) = 


d 

dt 


f 

J/=o 


Po v„ 


+ 


p cos 6 


L r (1 — M r ) cr(l— M r ) \ 


dS + 


ret 


f 

Jf=0 


p cos 6 

L r 2 "7T^A^") J 


dS 


(39) 


ret 


The observer time differentiation of Formulation 1 is always taken numerically in the applications of 
this result. 

Remark 5- The solution of the wave equation with derivatives acting on the inhomoge- 
neous source term 


Equations (5) and (6) are of the following types: 

j dq j dq , o N 

□ 0 = — and □ (pi = xe !R , q differentiable in x and t) 

dt dxi 

Let us consider the first equation. Assume that we find the solution of the following wave equation: 

□ 2 ip = q 

Then, by taking the time derivatives of both sides of eq. (41), we have the following result: 


d 


df dq 


— nr 0 = □ 

dt dt dt 

Comparing this equation with the first wave equation in eq. (40), we have 

df (x, t) 


(p (x, t) = 


dt 


(40) 


(41) 


(42) 


( 43 ) 


This is the result we used in writing down eq. (32). Similarly, we have 
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d 0 9 dib da 

— n 2 if/ = n 2 — - = — 

dxj dx{ dxj 

When this result is compared with the second wave equation in eq. (40), we conclude that 
dif/ (jc, t ) 


0i O, 0 = 


dxi 


(44) 


(45) 


This is the result used in eq. (34). 

We make two further important comments here. First deriving the results in eqs. (43) and (45) using the 
Green’s function solution of the wave equations of eq. (40) in the integral form as it is often done in 
books and technical journals, is difficult and involves much more algebraic manipulations. Second, if 
we relax the differentiability of the source function q(x, t ) in eq. (40), then the results expressed by eqs. 
(43) and (45) are valid if all the derivatives d / dt and d / dx t on the right sides of eqs. (40), (43), and 
(45) are treated as generalized derivatives (See References 5 and 6). In fact, it is always more conve- 
nient to work in the space of generalized functions because of their powerful operational properties 
(ease of exchanging limit processes without worrying about differentiability, etc.). As a matter of fact, in 
many practical problems, differentiability or even the continuity of the source term cannot be guaran- 
teed. However, using generalized functions can overcome most, if not all, the mathematical ambiguities 
and difficulties. (End of Remark 5) 

Remark 6- In eq. (37), we have taken the observer time derivative out of the integral sign after the 
second equality sign. But we note that the integral sign in this equation stands for a four dimensional 
integral as given by eq. (14). The upper limit of the integration with respect to the source time is t. A 
keen reader would recognize that the Leibniz rule of differentiation under an integral sign (see Refer- 
ence 8 below) must be used to establish the validity of the operation in eq. (37). We will now demon- 
strate how to do this. 


Let us use the notation Q - p n, r , . We start with the following integral appearing on the right of eq. 
(37): 


t d Cnt urn j a d 

I = — Q(y, T)S(f) — — dydr = 


dt 


f'f 

*7—oo *7 [R 


r dt 

Now take the time derivative inside the integral applying Leibniz rule: 

G O, t) 

'R 


GO, t) 6(f) ----- dydr 

oo *7 IR-^ r 


-f'f 

*7 — oo Jr* 


GO, t) d 5 (g) 

— 7 — dydr + 

r dt 


f 

*7ir 


6[f(y, t)\ 6 \-\dy 


(46) 


(47) 


We will now show that the second integral is zero. Lirst note that 6(r/c) = c6(r). So we must show that 
the following integral vanishes: 

r q o, o 

h O) = I d[f(y, t)] 5(r)dy (48) 

Jr 3 r 


We use the spherical polar coordinates (r, ip, d) with center at the observer. We have 
d y = r 2 sin d d r d p dd. Therefore, the above integral can be written as: 
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n 2 it /->oo 

I Q (y, t) 6[f(y, t)]r sind d(r - s) dr d<p dd = 0 (49) 

Jo 


The reason is that r 6 (r - s) - s 6 (r - e) and, thus, the above triple integral is of the order of s and goes 
to zero as e -> 0 . Therefore, we have shown that 


— 
at J 

foiy, T)S(f) - r dy dr = 

r Q O’, T) S(f ) E 

n 

dy dT 


(50) 


Note that even though r is not a function of the observer time, we like to always associate r with 6(g) 
and write the observer time derivative operating on 6(g) / r in the integrand of eq. (50). 

For Leibniz rule of differentiation under an integral sign, I recommend volume 1 of the following book: 

8- R. Courant: Differential and Integral Calculus, 2 volumes, Interscience Publishers, 1936. (Pub- 
lished also more recently in Wiley Classic Library Series) 

This book is considered by some mathematicians as perhaps the best calculus book of the twentieth 
century. A somewhat modernized version of this book is by R. Courant and F. John (in two volumes). 
Both versions are treasure-troves of the most useful calculus results for applications. Volume 2 of this 
book also includes all the fine points of the subjects of transformation of variables, n-dimensional vol- 
umes and surfaces that one needs to understand the present work. These are difficult subjects that are no 
longer emphasized in calculus courses. They are discussed by Courant without the use of excessive 
abstract language employed by some modern authors of books and articles on mathematics. One must 
pay careful attention to all of the footnotes in the above calculus book because of the many important 
examples and comments there. (End of Remark 6) 


3.4- Derivation of Formulation 1A 


Derivation of this result is based on Formulation 1. As will be seen, the discussions in Subsection 3.3 
make the derivation of this formulation an exercise in partial differentiation, albeit a fairly complicated 
one algebraically. 

Let us look again at eq. (30): 


E (x, t\ //) = 


Q C V, T) 


Lr(l- M r ) 


J ret 


Q (y, t) 
r (1 - M r ) 


(x; t], r) 


J T = T e 


(51) 


We have shown that only the emission time T e is a function of the observer time. Therefore, if we need 
to take the observer time derivative of the above expression, we must use the chain rule as follows: 


( 2-1 - 

8t (x, t\ if) d 

\dt) x 

dt dr _ 


(52) 


where r (x, t\ t]) is simply the solution of r - t + | x - y (ij, t)| / c = 0 . 

Remark 7- Whenever you see a partial derivative with respect to a variable, stop and ask yourself what 
other variables are kept fixed. The notation of partial differentiation can be very confusing when work- 
ing with the wave equation with sources on a moving surface. We note that: 
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1- In deriving this equation, we assume that the moving surface / = 0 is defined in a Lagrangian frame 
i] fixed to the surface. This frame is also called the blade-fixed frame. We have y = y (//, r), so that 
r = r(x, y (//, r)). Also note that the surface integration is most conveniently carried out in the //-frame, 
i.e., we are really integrating over the surface f (rj) - 0. 

2- We have also shown that the emission time is analytically written as follows r = r e (x, t; rj) although 
we cannot find this function explicitly for rotating blades. 

Therefore, in partial differentiation of the left side of eq. (52), we are keeping the variables (x, //) fixed. 
All these also mean that we need to find dr e / dt on the right side of eq. (52) keeping the variables (x, //) 
fixed. (End of Remark 7) 

Now we find dr / dt. We have shown above that the emission time satisfies the following equation: 

g = T-t + r(x, y(ij,T))/c = 0 (53) 


where // is a given fixed point on the moving surface. By taking partial derivative with respect to t of 
eq. (53), we get 


V dt / ( x,Tj ) 


1 (dr 
1 + - — 
c V dt 


(■ x , if) 


d T e 
dt 


1 - AT, 


(■ X , Jl) 


(£) =0 

V dt J(x,rj) 


(54) 


where M r is the Mach number of the point ij in the radiation direction at the time r. Here we have used 
the following relations: 


dr 

dt 

dr 

dr 


( X , T]) 


dr ( dy { N 


dr 

dt 


(X, Tj) 


(x, rf) 


dy t \ dr 


= - n v: = - v,. 


(55) 


(56) 


where r, is the component of unit radiation vector and vy is the velocity of the point // in the radiation 
direction. 


From eqs. (54-56), we find 



(57) 


This result leads to the following important formula: 


d d [ 1 dq{ x, v, r) 

— [q(x, y, r)] = — [q(x, y (tj, t), i e (x, t; //))] = — 

dt m dt LI -M r dr 


(58) 


J ret 


We warn the readers once more that the right side of eq. (58) is a shorthand notation which must be 
interpreted as: 
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1 

dq{x , y, r) 


1 

dq(x, y (tj, t), r) 

A -M r 

dr 

ret 

1 - M r (x\ TJ, r) 

dT 


(59) 


This analytic expression seems complicated and intimidating. But the physical interpretation is simple if 
we work in the ij -frame but continue using the notations we use in the y-frame to reduce confusion. 
Equation (59) can be simply rewritten as: 


1 

dq(x, y, r) 


1 

dq(x, tj, r) " 

A -M r 

dz 

ret 

1 - M r (x\ TJ, r) 

dz 


Since we perform noise prediction by dividing a blade surface into panels and compute the contribu- 
tions of the sources on each panel separately, we are indeed working with a fixed t] for each panel. 
Therefore, the physical meaning of every term in eq.(60) is quite clear and pertain to physical properties 
of a given panel at the emission time. In particular, when we deal with the loading source of Formula- 
tion 1A, p (//. r) is precisely the unsteady blade surface gage pressure that a transducer fixed to the 
blade surface at position rj measures. Therefore, in this formulation p = d p (rj, t) / dr . We discuss these 
ambiguities in notation because in some of our other formulations, the surface pressure p stands for 
p(y, t) which is the unsteady gage pressure measured by a transduced fixed to the undisturbed medium. 
Note that p (tj, z) = p(y (tj, r), r) and, if we had not abused the notation, we should have used 
pin, t) = p(y(n , t), r). 


We are now just one short step away from Formulation 1A. We use Formulation 1, eq. (39), and take the 
observer time derivative inside the first integral to obtain: 


4 71 p' (x, t) = 


Jf= 


d 


Po A 


+ 


p cos 6 


f=0 dt Lr(l - M r ) cr (1 - M r ) \ 


dS + 


ret 


f [ 

J/=o L 


pcosO 
r 2 (1 - M r ) _ 


dS 


ret 


(61) 


Next we use eq.(60) inside the first integral to get 

4/r p' (x, t ) = 4 n{p' T (x, t ) + p' L (x, t)) = 


U 


1 


d 


Po v r 


+ 


p cos 6 


1 -M r dr L r (1 - M r ) cr(l-M r )\ 
dS 


dS + 


f 

p cos 6 

)f= 0 

_ r 2 (1 - M r ) _ 


(62) 


Now we must remember that everything inside the integrand of the first integral, before the variable t 
is replaced with r e , is a function of either (//. r) or (x; //. r) . First note that cos 6 = n- r and M, - M • r . 
The variables v n , M, p, andn are functions of (//. r), and the variables randr are functions of (x; //. r). 
To evaluate the integrand of the first integral of eq. (62), we need the source time derivatives of all these 
seven variables as we will see below. 

We have already talked about the meaning of the time derivative of the surface pressure p .The source 
time derivative of v n will be derived next. We have the following algebra 
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V„ = 


d 

dr 


(v • n) = v-n + V'h = a n + v • ajxn 


(63) 


where a = v-dv/dr-cM is the acceleration of the point ij and oj is the angular velocity of the blade 
surface. Note that, in general, v-cjxn + 0 on a blade. We can now easily show the following results: 


dr _ dr dy { (tj, r) 
dr dyt dr 


(64) 


dfj _ fj v r - Vj 
dr r 


(65) 


dM r 

dr 


1 

c r 


(r f v/ + v 2 - v 2 ) = fi Mi + 


c(M 2 -M 2 ) 

r 


( 66 ) 


We can next perform the tedious source time differentiation in eq.(62) and utilize the above results in 
the algebra. Here we separate the thickness and loading noise results as in Reference 4. The Thick- 
ness and Loading noise components of Formulation 1A is 


47r p' T (x, t) = 

J' \ Po v n 


+ 


Po v n r, Mi 


/=o L r (1 - M r ) r (1 — M r )~ 


dS + 


Jret 


f 

Jf=0 


Po c V n (M r -M 2 ) 
r 2 (1 -M r ) 3 


dS 


(67) 


Jret 


p' L (x,t)= I 

Jf= I 

f 

J/=o 


p cos 6 


+ 


Mi p cos 6 


dS + 


f=o[cr(l-M r ) cr{\-M r )~ 
p (cos 6 - Mj rij) (My - A/ 2 ) p cos 9 


Jret 


r 2 (1 -M r Y 


r 2 (1 — M,-Y 


dS 


Jret 


( 68 ) 


We have written these equations differently from those in Reference 4 by separating the near field 
terms (order l/r 2 ) from the far field terms (order 1/r). We have left out much algebraic manipula- 
tions in obtaining the above results. The derivation of Formulation 1A has been checked many times by 
acoustic researchers. We mention here that although the combination of variables My - M 2 appears in 
eq.(66), the appearance of the combination of variables M r - M 2 in eqs.(67) and (68) is correct. 

Remark 8- The loading noise part of Formulation 1A can be derived in other ways. For example, we 
have 
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4 71 p' L (X, t) = 


8 


7 


p cos 6 


-L 


8 


dxi J f=0 L r( 1 - Af r ) 
p cos 0 


(x; //, r) 


rf5 = 


=o dx f L r (1 - M r ) 


(x; T], r) 


>lf-0 ^ dXi 


p cos 6 
L riX-Mr) 


J r = r e (x, ?;»/) 


(x; T], r) 


J r = r f (x, t;ij) 

dS = 

+ 


dr ( d \ pcosO 


dxj V dr / (JC;J|) L r (1 — M r ) 


(x; //, r) 




T = T C (X, t\7]) 


Now, From eq. (53), we find 


dr 

d Xj 


+ — - M r 


8t 

dxj 


= 0 


From this, we find 

dr n 

dxi c(l—M r ) 


( 69 ) 


( 70 ) 


( 71 ) 


Use this result in eq. (69) and carry out all the differentiation with respect to Xj to get the loading part of 
Formulation 1A. Note that in this derivation, we do not need to use the Leibniz rule of differentiation 
under the integral sign when we take the observer space derivative inside of the integral in eq. (69). 

We did not follow this procedure originally because we already had Formulation 1 in our possession and 
only the observer time differentiation had to be performed analytically for both thickness and loading 
noise results. There were some unexpected applications of the procedure of converting the observer 
space derivative to observer time derivative. Farassat and Brentner used this procedure to derive a 
formulation for quadrupole noise prediction similar to Formulations 1 which should properly have been 
called Formulation Q1 (Reference 9, eq. (14), Reference 10, eq. (4)). Later, Brentner derived a sec- 
ond formulation for quadrupole noise prediction similar to Formulation 1A called Formulation Q1A 
(Reference 10, eq. (17)) . This latter result has been successfully used in helicopter rotor noise predic- 
tion. The primary references to these works are: 

9- F. Farassat and Kenneth S. Brentner: The uses and abuses of the acoustic analogy in helicopter 
rotor noise prediction, Journal of the American Helicopter Society, 1988, 33, 29-36 

10- K. S. Brentner:An efficient and robust method for predicting helicopter high-speed impulsive 
noise, Journal of Sound and Vibration, 1997, 203(1), 87-100 

Had we not developed the analytic technique of converting the observer space derivative exactly to the 
observer time derivative to derive the loading noise component of Formulation 1, we probably would 
not have been able to obtain our quadrupole noise formulations. (End of Remark 8) 
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3.5- Moving Observer 

In most applications, we want to compute the noise in a frame moving with the aircraft, i.e., as measured 
by a microphone attached to the aircraft. All propeller and helicopter rotor noise prediction codes at 
Langley Research Center have this capability. It is important to realize that the observer variable x in 
Formulations 1 and 1 A is in the frame fixed to the medium. Let us attach a frame fixed to the aircraft 
with axes parallel to the medium-fixed x- frame. We will call this A-frame. Assume that at the time t = 0 
the two frames coincide. Let the velocity of the aircraft be V(t). Then the origin of the A-frame will be 
at the point 

X 0 (t)= f V(t')dt' (72) 

Jo 


Now a point A in the A-frame will be at the point x = A + AqL) at the time t. Therefore, if we want to 
find the acoustic pressure p' (A, t) as a function of time in the A-frame from our formulations, we use 
the following formula: 


~P (*, t) = p' (X + X 0 (t), t) 


(73) 


The interpretation of this result is simple. To calculate the acoustic pressure p (A, t ) in the aircraft fixed 
frame, make sure you compute first where the observer (or the microphone) is in the x-frame at the time 
t. Then use that observer position in the formulations given by eq. (73). 


4- Concluding Remarks 

Langley Research Center researchers have been involved in the prediction of the noise of rotating blades 
since the nineteen thirties. It was after the introduction of high speed digital computers in nineteen 
sixties that it was possible to use realistic geometry and kinematics of the rotors and propellers in the 
noise prediction process. The pace of aeroacoustic research increased since the early nineteen seventies 
because of the public pressure to reduce aircraft noise particularly around airports. Both the govern- 
ments and the aircraft industry around the world realized that aeroacoustic research can lead to substan- 
tial aircraft noise reduction. The invention of high speed digital computers played a major role in this 
effort in the areas of aeroacoustic modeling, transducer design, data analysis, full scale flight and wind 
tunnel tests. Advanced mathematics have played a vital role in all these areas although this fact is not 
mentioned often by the researchers. 

In developing noise prediction models at Langley Research Center for propellers and helicopter rotors, 
we used primarily the acoustic analogy based the FW-H equation. We developed general solutions of 
this equation which allowed using realistic blade geometry and kinematics. We also did not want to be 
limited to problems that generated periodic acoustic waveforms. Finally, we required that our formula- 
tions be valid in both near and far fields. These conditions could be met if we worked in the time 
domain and then applied a Fourier transform in time to obtain frequency domain results. This approach 
has been very fruitful so far and is being followed by many other researchers around the world. 

The current trend in aircraft noise prediction appears to be toward using FW-H pc j s for all helicopter rotor 
and propeller noise prediction. This would make Formulation 1A even more useful in the future for all 
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rotating blade noise prediction problems. For example, it makes sense to use this formulation for predict- 
ing supersonic propeller noise when an advanced unsteady aerodynamic code becomes available. This is 
what we are proposing at present for this problem. 

Now I would like to give some practical advice for developing a noise prediction code based on my own 
experience at Langley Research Center: 

1- Always derive fully the acoustic formulation you want to use. This way you understand the exact 
meaning of all the terms in the formulation as well as its subtleties, e.g., the fluid velocities are defined 
relative to what frame. Also, it is possible that some misprints in the printed reference, such as a sign 
error, can cause serious errors in acoustic calculations. 

2- Spend a lot of time designing the algorithms you want to use and studying the flow of information in 
the code. This is the most important step in the reduction of computation time and it should be given 
serious attention. We give some examples below: 

-Use advanced surface integration techniques such as Gauss-Legendre integration . 

- Emission time computation can be very time consuming on a computer. A well-designed algorithm 
here is essential for an efficient code. 

- Sometimes a do loop in the wrong place in the code can increase the computation time considerably. 
This is where the experience of the code developer becomes important. 

3- Since we have derived exact results here, the readers should realize that in the near field the Fourier 
transform of the computed acoustic pressure can have a large constant term. This term is negative on the 
suction side of the rotor or propeller disc and positive on the pressure side as the physics of the problem 
would dictate. This is confusing to a novice in the field of aeroacoustics because the measured results 
generally do not show this constant value (called DC shift by experimenters). Most commonly used 
microphones cannot measure this shift. These microphones are said to be AC coupled and remove the 
DC shift. 
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